import re  # regular expression 正则表达式
import numpy as np


class Geometry:
    """
    抽象几何类型
    """

    def __init__(self):
        self.dimension = None  # 几何对象的维数
        self.bbox = None

    def geometryType(self):
        """
        返回该对象的几何类型（str）

        :return: str
        """
        pass

    def asText(self):
        """
        将几何对象转为wkt文本描述（str）

        :return:
        """
        pass

    def fromText(self, wkt, tagged=True):
        """
        解析wkt文本得到几何对象（Geometry）

        :param wkt:
        :param tagged: 是否有类型标记("point" "linestring" "polygon")
        :return:
        """
        pass

    def fromXY(self, *args):
        """
        输入坐标构造几何类型。

        :param args: 参数元组(tuple)
        """
        pass

    def envelope(self):
        """
        返回该几何体的最小外包矩形，是一个几何对象（geometry）

        :return:
        """
        pass


class Point(Geometry):  # 继承自()内的类
    def __init__(self):
        super().__init__()
        self.dimension = 2  # 点是二维的
        self.x = None  # 点的x和y
        self.y = None

    def fromText(self, wkt, tagged=True):  # tagged参数默认为True，为false时，跳过if判断
        if tagged and wkt[:5].lower() != "point":  # wkt前5个字符为point, 可以将常量("point")写在等号左边
            return
        # 用正则表达式匹配所有数字
        # main里的point的wkt例子 'point(-1 3.14)'
        # ? 前面的字符最多出现一次（即0或1次）, -?代表x或y的坐标值可以是正的，也可以是负的
        # + 前面的字符至少出现一次
        # \d 匹配数字
        # \.\d 匹配得到小数点后的数字
        # 返回的nums类型为列表
        nums = re.findall(r'-?\d+(?:\.\d*)?', wkt)
        self.x = float(nums[0])  # x和y的类型变为float
        self.y = float(nums[1])
        self.bbox = [self.x, self.x, self.y, self.y]

    def asText(self):
        return f"Point({self.x} {self.y})"  # 加f后可以在字符串里使用花括号括起来的变量和表达式

    def fromXY(self, x, y):
        self.x = x
        self.y = y
        self.bbox = [self.x, self.x, self.y, self.y]

    def geometryType(self):
        return "Point"

    def envelope(self):  # 点的envelope是点。应该不需要用到点的envelope
        return self

    def equals(self, geom):  # 是否和另一个几何对象完全一样
        flag = False
        if geom.geometryType() == self.geometryType():
            if geom.X() == self.x and geom.Y() == self.y:
                flag = True
        return flag

    def X(self):
        return self.x

    def Y(self):
        return self.y

    def distance_to_line(self, p1, p2):  # 点到线段p1p2的距离
        # 注意 linestring.points[i]得到的是一个含x和y的tuple元组（大小为2）
        # 即 p1 是一个这样的元组 (x, y)
        # 使用numpy进行计算
        # 先定义三个向量
        p1p2 = np.array([p2[0] - p1[0],
                        p2[1] - p1[1]])
        p2p1 = -p1p2
        p1self = np.array([self.x - p1[0],
                          self.y - p1[1]])

        p2self = np.array([self.x - p2[0],
                          self.y - p2[1]])

        # 向量p1p2和p1self的夹角 = 向量点乘 / 向量的模
        costheta = np.dot(p1p2, p1self) / \
                (np.linalg.norm(p1p2) * np.linalg.norm(p1self))
        thetap1 = np.arccos(costheta)
        dist_vertical = np.linalg.norm(p1self) * np.sin(thetap1) # 点到线段的垂直距离
        print("第一个夹角为 " + str(np.rad2deg(thetap1)))

        costheta = np.dot(p2p1, p2self) / \
                (np.linalg.norm(p2p1) * np.linalg.norm(p2self))
        thetap2 = np.arccos(costheta)
        print("第二个夹角为 " + str(np.rad2deg(thetap2)))

        # arr = np.array([[self.x, self.y],
        #                [p1[0], p1[1]],
        #                [p2[0], p2[1]]])

        # 如果两个角都是锐角，则取垂直距离
        if np.rad2deg(thetap1) <= 90 and np.rad2deg(thetap2) <= 90:
            dist = dist_vertical
        else:  # 否则取点到两个端点的距离最小值
            dist = min(np.linalg.norm(p1self),
                       np.linalg.norm(p2self))
        print("距离为 " + str(dist))
        return dist

    def distance(self, geom):
        if geom.geometryType() == "Point":
            dist = np.sqrt((self.x - geom.X()) ** 2 +
                           (self.y - geom.Y()) ** 2)
            return dist
        if geom.geometryType() == "LineString":
            # 距离 = min(点到每一条线段的距离)
            dist = []  # 距离列表
            for i in range(len(geom.points) - 1):  # 3个点构成两条线段 -1
                temp_dist = self.distance_to_line(geom.points[i],
                                                  geom.points[i+1])
                dist.append(temp_dist)
                print("点到第" + str(i+1) + "条线段的距离为 " + str(temp_dist))
            return min(dist)  # 返回距离列表的最小值

    def intersects_two_line(self, p2, p3, p4):  # 判断线段self p2  与 线段p3p4是否相交
        arr_point = np.array([self.x, self.y])
        A = self
        B = p2
        C = p3
        D = p4
        # 判断线段AB  与 线段CD是否相交
        AC = np.array([C.x - A.x,
                       C.y - A.y])
        AD = np.array([D.x - A.x,
                       D.y - A.y])
        BC = np.array([C.x - B.x,
                       C.y - B.y])
        BD = np.array([D.x - B.x,
                       D.y - B.y])
        CA = -AC
        CB = -BC
        DA = -AD
        DB = -BD

        # 向量叉积x 的乘积 异号
        # 这个算法在参考博客的最后
        if np.cross(AC, AD) * np.cross(BC, BD) <= 0 \
                and np.cross(CA, CB) * np.cross(DA, DB) <= 0:
            return True
        else:
            return False

    def intersects_linearRing(self, l):  # 射线法: 返回点和线环的交点数量
        x_max = np.inf  # x轴的正无穷
        # https://blog.csdn.net/daming98/article/details/79561777
        num_points = 0  # 交点数
        arr_point = np.array([self.x, self.y])
        for i in range(len(l.points) - 1):  # 3个点构成两条线段 -1
            p1 = l.points[i]
            p2 = l.points[i+1]
            # 判断(self.x, self.y)和(+oo, self.y)组成的线段AB  与 线段p1p2(记为CD)是否相交
            AC = np.array([p1[0] - self.x,
                          p1[1] - self.y])
            AD = np.array([p2[0] - self.x,
                          p2[1] - self.y])
            BC = np.array([p1[0] - x_max,
                          p1[1] - self.y])
            BD = np.array([p2[0] - x_max,
                          p2[1] - self.y])
            CA = -AC
            CB = -BC
            DA = -AD
            DB = -BD

            # 向量叉积x 的乘积 异号
            # 这个算法在参考博客的最后
            if np.cross(AC, AD) * np.cross(BC, BD) <= 0 \
                    and np.cross(CA, CB) * np.cross(DA, DB) <= 0:
                num_points += 1
            print(num_points)
        return num_points

    def intersects(self, geom, tolerance=0):
        """
        判断自身与geom是否在给定容差下相交，若相交返回True，否则返回False

        :param geom: 输入的几何图形
        :param tolerance: 容差
        :return: bool
        """
        # 点与点或线相交。考虑tolerance，即距离小于tolerance认为是相交。
        #  距离计算建议用矢量法。https://www.cnblogs.com/lyggqm/p/4651979.html
        if geom.geometryType() == "Point" or geom.geometryType() == "LineString":  # 点与点或线相交
            distance = self.distance(geom)
            if distance <= tolerance:
                print("相交")
                return True  # 考虑tolerance，即距离小于tolerance认为是相交
            else:
                print("不相交")
                return False
        # 点与多边形相交。不考虑tolerance
        # 计算点的射线（y不变，x向正无穷大）与外环和所有内环的交点和
        elif geom.geometryType() == "Polygon":
            num_point = 0
            num_point += self.intersects_linearRing(geom.exteriorRing())
            for i in range(geom.numInteriorRing()):  # 遍历内环
                num_point += self.intersects_linearRing(geom.interiorRingN(i))

            # 如果点的坐标含小数，不是整数，射线法就不用判断特殊情况了
            if 0 == num_point % 2:  # 交点为偶数，不在多边形内
                print("不相交")
                return False
            else:  #交点为奇数，在多边形内
                print("相交")
                return True

        else:
            print("不相交")
            return False


class Curve(Geometry):
    def __init__(self):
        super().__init__()
        self.dimension = 2
        # 声明待赋值的字段：点的列表，bbox和4个最大最小值
        self.points = None
        self._xMax = None
        self._xMin = None
        self._yMax = None
        self._yMin = None
        self.bbox = None

    def length(self):
        pass

    def startPoint(self) -> Point:
        pt = Point()
        pt.fromXY(self.points[0][0], self.points[0][1])
        return pt

    def endPoint(self) -> Point:
        pt = Point()
        pt.fromXY(self.points[-1][0], self.points[-1][1])
        return pt

    def isClosed(self) -> bool:
        firstPt = self.startPoint()
        lastPt = self.endPoint()
        if firstPt.equals(lastPt):
            return True
        else:
            return False


class LineString(Curve):
    def fromText(self, wkt, tagged=True):
        # 参考：https://stackoverflow.com/questions/16731461/parsing-a-wkt-file
        # tagged为true 并且 wkt类型不为LineString 则 return
        if tagged and wkt[:10].lower() != "linestring":
            return
        nums = re.findall(r'-?\d+(?:\.\d*)?', wkt)  # 用正则表达式匹配所有数字
        self.points = []  # 坐标列表

        self.points = list(zip(*[iter(map(float, nums))] * 2))  # 坐标列表。格式：[(x1,y1), (x2,y2)...]

        # 求坐标的范围。numpy的计算效率远高于for循环。
        coords_arr = np.array(self.points)
        self._xMax = np.max(coords_arr[:, 0])  # 第0列的最大值
        self._xMin = np.min(coords_arr[:, 0])
        self._yMax = np.max(coords_arr[:, 1])
        self._yMin = np.min(coords_arr[:, 1])
        self.bbox = [self._xMin, self._xMax, self._yMin, self._yMax]

    def length(self):
        p = self.points
        len_list = [np.sqrt((p[i][0] - p[i - 1][0]) ** 2 +
                            (p[i][1] - p[i - 1][1]) ** 2) for i in range(1, len(p))]
        sum_len = sum(len_list)  # 求和计算曲线长度
        return sum_len

    def asText(self):
        # 一行代码太长时，敲enter，自动生成\(代表换行)
        # LineString main中的wkt示例 'LINESTRING(-1.1 0, 1 1, 3 -1)'
        return "LineString(" + \
               ",".join([f"{xy[0]} {xy[1]}" for xy in self.points]) + \
               ")"

    def fromXY(self, coords):
        self.points = coords
        coords_arr = np.array(coords)
        # 求坐标的范围
        self._xMax = np.max(coords_arr[:, 0])
        self._xMin = np.min(coords_arr[:, 0])
        self._yMax = np.max(coords_arr[:, 1])
        self._yMin = np.min(coords_arr[:, 1])
        self.bbox = [self._xMin, self._xMax, self._yMin, self._yMax]

    def envelope(self):
        bbox = Polygon()
        bbox.fromXY([(self._xMin, self._yMin),
                     (self._xMin, self._yMax),
                     (self._xMax, self._yMax),
                     (self._xMax, self._yMin),
                     (self._xMin, self._yMin)])
        return bbox

    def geometryType(self):
        return "LineString"

    def numPoints(self):
        """
        返回线段上点的数量
        :return:
        """
        return len(self.points)

    def pointN(self, n) -> Point:
        """
        返回点列表中第n个Point。 注意：列表中n从0开始

        :param n:
        :return:
        """
        pt = Point()
        pt.fromXY(self.points[n][0], self.points[n][1])
        return pt

    def intersects(self, geom, tolerance=0):
        if geom.geometryType() == "Point":
            return geom.intersects(self, tolerance)
        elif geom.geometryType() == "LineString":  # 线与线相交
            # 暂时没考虑包围盒快速过滤操作
            # 分别对两个linestring的所有线段进行判断
            for i in range(self.numPoints() - 1):
                for j in range(geom.numPoints() - 1):
                    p1 = self.pointN(i)
                    p2 = self.pointN(i+1)
                    p3 = geom.pointN(j)
                    p4 = geom.pointN(j+1)
                    if p1.intersects_two_line(p2, p3, p4):
                        print("线与线相交")
                        return True
            print("线与线不相交")
            return False
        # 线与多边形相交
        elif geom.geometryType() == "Polygon":  # LineString.intersects(Polygon)
            # 1. 线与geom多边形的外环是否相交
            if self.intersects(geom.exteriorRing()):
                print("线与多边形相交")
                return True
            else:
                # 是否与多边形的内环相交
                for i in range(geom.numInteriorRing()):
                    if self.intersects(geom.interiorRingN(i)):
                        return True
                # 2. 点是否在多边形内部（只用判断一个点即可）
                if self.pointN(0).intersects(geom):
                    for i in range(geom.numInteriorRing()):
                        if self.pointN(0).intersects(geom.interiorRingN(i)):  # 在内环以内
                            return False
                    return True
                else:  # 3. 不符合1,2 则不相交
                    return False
        else:
            return False


class LinearRing(LineString):
    def fromText(self, wkt, tagged=True):
        super().fromText(wkt, tagged)
        assert self.isClosed(), "环未闭合"

    def fromXY(self, coords):
        super().fromXY(coords)
        assert self.isClosed(), "环未闭合"


class Surface(Geometry):
    pass


class Polygon(Surface):
    def __init__(self):
        super().__init__()
        self.dimension = 2
        # 声明待赋值的字段：外环、内环和bbox
        self.__exteriorRing = None
        self.__interiorRings = None
        self.bbox = None

    def bbox(self):
        return self.__exteriorRing.bbox()

    def geometryType(self):
        return "Polygon"

    def fromText(self, wkt, tagged=True):
        if tagged and wkt[:7].lower() != "polygon":
            return
        ringStrList = re.findall(r'\([^()]+\)', wkt)
        ringList = []
        for ringStr in ringStrList:
            linearRing = LinearRing()
            # tagged=False，使不含linestring的wkt能被解析得到几何对象
            linearRing.fromText(ringStr, tagged=False)
            ringList.append(linearRing)
        # 外环
        self.__exteriorRing = ringList[0]
        # 内环
        self.__interiorRings = ringList[1:]
        self.bbox = self.__exteriorRing.bbox

    def asText(self):
        # 外环  用[10:]删除wkt里的linestring这10个字符
        exteriorStr = self.__exteriorRing.asText()[10:]
        # 内环
        interiorStr = ""
        for interior in self.__interiorRings:
            interiorStr += "," + interior.asText()[10:]
        return f"Polygon({exteriorStr}{interiorStr})"

    def fromXY(self, coords):
        # 暂时只考虑没有洞的多边形
        self.__exteriorRing = LinearRing()
        self.__exteriorRing.fromXY(coords)
        self.__interiorRings = []
        self.bbox = self.__exteriorRing.bbox

    def envelope(self):
        return self.__exteriorRing.envelope()

    def exteriorRing(self):
        """
        获取外环

        :return:
        """
        return self.__exteriorRing

    def numInteriorRing(self):
        """
        获取内环个数

        :return:
        """
        return len(self.__interiorRings)

    def interiorRingN(self, n):
        """
        获取下标为n的内环。内环下标范围：[0, numInteriorRing()-1]

        :param n: 要获取的内环下标
        :return:
        """
        assert n < len(self.__interiorRings)
        return self.__interiorRings[n]

    def intersects(self, geom):
        if geom.geometryType() == "Point":
            return geom.intersects(self, 0)
        elif geom.geometryType() == "LineString":
            return geom.intersects(self, 0)
        # 多边形与多边形相交
        elif geom.geometryType() == "Polygon":
            # 注：此处暂时只考虑了self没有内环的情况
            # 1. 先判断线与线有没有相交，有就返回True；
            # 2. 继续判断self的某一点在不在框选矩形内部，在就返回True
            # 3. 反过来判断geom的某一个点是否在self内部，在就返回True geom在内环内不行
            # 4. 返回False
            if self.intersects(geom.exteriorRing()):
                print("多边形与多边形相交")
                return True
            else:
                # 判断是否与内环相交
                for i in range(geom.numInteriorRing()):
                    if self.intersects(geom.interiorRingN(i)):
                        return True

                # 判断self的点在geom内部还是外部
                pt = self.exteriorRing().pointN(0)
                exteriorPolygon = Polygon()
                exteriorPolygon.fromXY(geom.exteriorRing().points)
                if pt.intersects(exteriorPolygon):
                    flag = True
                    for i in range(geom.numInteriorRing()):
                        interiorPolygon = Polygon()
                        interiorPolygon.fromXY(geom.interiorRingN(i).points)
                        if pt.intersects(interiorPolygon):  # 在内环内部，说明不相交
                            flag = False
                            break
                    return flag

                else:  # 在外部
                    return False
        else:
            return False


def test_intersect_point():
    pt_wkt = "point(1 3)"
    pt = Point()
    pt.fromText(pt_wkt)
    print("点为" + pt_wkt)

    linestring_wkt = 'LINESTRING(-1.1 0, 1 1, 3 -1, 4 6)'
    ls = LineString()
    ls.fromText(linestring_wkt)
    print("线为" + linestring_wkt)

    polygon_wkt = 'POLYGON((0 0,0.5 4, 4 4, 3.5 -1, 0 0),(1 1, 2 1, 2 2, 1 2, 1 1))'
    polygon = Polygon()
    polygon.fromText(polygon_wkt)
    print("多边形为" + polygon_wkt)

    #result = pt.intersects(ls, 1)
    #result = pt.intersects(ls, 3)
    pt.intersects(polygon)

    pt.fromXY(1.5, 1.5)
    print("点为" + str(pt.x) + " " + str(pt.y))
    pt.intersects(polygon)

    pt.fromXY(-1, 1.5)
    print("点为" + str(pt.x) + " " + str(pt.y))
    pt.intersects(polygon)

def test_intersect_linestring():
    linestring_wkt = 'LINESTRING(0 0, 2 2, 2 0)'
    ls1 = LineString()
    ls1.fromText(linestring_wkt)
    print("线1为" + linestring_wkt)

    linestring_wkt = 'LINESTRING(0 1, 2 3, 3 3)'
    ls2 = LineString()
    ls2.fromText(linestring_wkt)
    print("线2为" + linestring_wkt)

    polygon_wkt = 'POLYGON((-1 -1, 3 -1, 3 3, -1 3, -1 -1))'
    polygon = Polygon()
    polygon.fromText(polygon_wkt)
    print("多边形为" + polygon_wkt)

    ls1.intersects(ls2)
    print('\n')
    ls1.intersects(polygon)
    print('\n')

    polygon_wkt = 'POLYGON((1 1, 3 1, 3 3, 1 3, 1 1))'
    polygon = Polygon()
    polygon.fromText(polygon_wkt)
    print("多边形为" + polygon_wkt)

    ls1.intersects(polygon)

def test_intersect_polygon():

    polygon_wkt = 'POLYGON((-1 -1, 3 -1, 3 3, -1 3, -1 -1))'
    polygon1 = Polygon()
    polygon1.fromText(polygon_wkt)
    print("多边形1为" + polygon_wkt)

    polygon_wkt = 'POLYGON((2 2, 4 2, 4 4, 2 4, 2 2))'
    polygon2 = Polygon()
    polygon2.fromText(polygon_wkt)
    print("多边形2为" + polygon_wkt)

    polygon1.intersects(polygon2)

    print('\n')

    polygon_wkt = 'POLYGON((4 4, 5 4, 5 5, 4 5, 4 4))'
    polygon3 = Polygon()
    polygon3.fromText(polygon_wkt)
    print("多边形3为" + polygon_wkt)

    polygon1.intersects(polygon3)

def test_distance_point_to_line():
    pt_wkt = "point(1 3)"
    pt = Point()
    pt.fromText(pt_wkt)
    print("点为" + pt_wkt)

    linestring_wkt = 'LINESTRING(-1.1 0, 1 1, 3 -1, 4 6)'
    ls = LineString()
    ls.fromText(linestring_wkt)
    print("线为" + linestring_wkt)

    print("点到折线的距离为" + str(pt.distance(ls)))

def test_geometry():
    # point
    pt_wkt = "point(-1 3.14)"
    pt = Point()
    pt.fromText(pt_wkt)
    print(pt.X(), pt.Y())
    pt.fromXY(3, 1.1)
    print(pt.X(), pt.Y())

    # linestring
    linestring_wkt = 'LINESTRING(-1.1 0, 1 1, 3 -1)'
    ls = LineString()
    ls.fromText(linestring_wkt)
    print(ls.points)
    print(ls.bbox)
    print(ls.length())
    print(ls.envelope().exteriorRing().points)  # envelope返回的是Polygon，不是bbox
    print(ls.startPoint().X(), ls.startPoint().Y())
    print(ls.isClosed())

    # polygon
    polygon_wkt = 'POLYGON((0 0,0.5 4, 4 4, 3.5 -1, 0 0),(1 1, 2 1, 2 2, 1 2, 1 1))'
    polygon = Polygon()
    # fromText
    polygon.fromText(polygon_wkt)
    print(polygon.bbox)
    # exteriorRing、points
    print("polygon外环的坐标：", polygon.exteriorRing().points)
    # envelope
    print("polygon的envelope的坐标：", polygon.envelope().exteriorRing().points)
    # numInteriorRing
    print("polygon的内环的数量：", polygon.numInteriorRing())
    # fromXY
    polygon.fromXY([(0, 0), (1, 3), (3, -1), (0, 0)])
    # asText
    print(polygon.asText())

    polygon_wkt = 'POLYGON((0 0,0 4, 4 4, 4 0))'
    try:
        polygon.fromText(polygon_wkt)  # 将报错：环未闭合
    except Exception as e:
        print(str(e))


if __name__ == "__main__":
    test_geometry()
    #test_distance_point_to_line()
    #test_intersect_point()
    #test_intersect_linestring()
    #test_intersect_polygon()
